Q1 and Q2 duration in stress and ictus

library("lme4")

Loading required package: Matrix

library("ggplot2")
library("broom")
library("dplyr")

Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

data_df <- read.csv("vowel_master.csv", header= TRUE)


data_df$ictus <- as.factor(data_df$ictus)
data_df$stressed <- as.factor(data_df$stressed)
data_df$stressed <- recode_factor(data_df$stressed, "0" = "unstressed","1" = "stressed")
data_df$quantity <- as.factor(data_df$quantity)
data_df$song <- as.factor(data_df$song)
data_df$word <- as.factor(data_df$word)
#head(data_df)

#to look at vowel duration in song-level and word-level prominence positions, I filter out only syllables of Q1 and Q2, and only those shapes which surface in both stressed and unstressed positions 
strictus_data <- filter(data_df, data_df$quantity!="3" & data_df$shape != "VCC" & data_df$shape != "VVC" & data_df$shape != "VV")

#strictus_data <- filter(strictus_data, nchar(strictus_data$segment)==1)
head(strictus_data)

NA
NA

ggplot(strictus_data,aes(ictus,duration)) + 
  geom_jitter(color = "violet", alpha = 0.5) +
  facet_wrap("song") +
  labs(x = "beat position", y = "mean vowel duration (ms)") + 
  theme_minimal()
ggsave("song_icK_dur.png")

Saving 7.06 x 4.36 in image

  
ggplot(strictus_data,aes(stressed,duration)) + 
  geom_jitter(color = "orange", alpha = 0.5) +
  facet_wrap("song") + 
  labs(x = "word stress", y = "within-song mean vowel duration (ms)") + 
  theme_minimal()
ggsave("song_str_dur.png")

Saving 7.06 x 4.36 in image

ggplot(strictus_data,aes(ictus,duration)) + 
  geom_jitter(color = "violet", alpha = 0.5) +
  facet_wrap("performer") + 
  labs(x = "beat position", y = "mean vowel duration (ms)") + 
  theme_minimal()
ggsave("perf_ick_dur.png")

Saving 7.06 x 4.36 in image

ggplot(strictus_data,aes(stressed,duration)) + 
  geom_jitter(color = "orange", alpha = 0.5) +
  facet_wrap("performer") + 
  labs(x = "word stress", y = "within-subject mean vowel duration (ms)") + 
  theme_minimal()
ggsave("perf_str_dur.png")

Saving 7.06 x 4.36 in image

interaction.plot(strictus_data$quantity,strictus_data$ictus,strictus_data$duration,xlab="quantity",ylab="mean vowel duration(ms)",trace.label = "beat")

interaction.plot(strictus_data$stressed,strictus_data$ictus,strictus_data$duration,xlab="stressed",ylab="mean vowel duration(ms)",trace.label = "beat")

ggplot(strictus_data,aes(x=stressed ,y=duration,color = ictus)) + 
  geom_violin() + stat_ydensity() + theme_minimal() +
  labs( xlab = "word stress" ,ylab = "vowel duration (ms)") + 
  facet_wrap("quantity", labeller = label_both)
  
ggsave("dur_density_qfac.png")

Saving 7.06 x 4.36 in image

ggplot(strictus_data,aes(x=shape,y=duration, color = ictus)) +
geom_violin() + stat_ydensity() + theme_minimal() +
labs( xlab = "syllable shape" ,ylab = "vowel duration (ms)") 

ggsave("dur_shapeictus.png")

Saving 7.06 x 4.36 in image

ggplot(strictus_data,aes(x=shape,y=duration, color = stressed)) +
geom_violin() + stat_ydensity() + theme_minimal() +
labs( xlab = "syllable shape" ,ylab = "vowel duration (ms)") 

ggsave("dur_shapestress.png")

Saving 7.06 x 4.36 in image

lindesignmod <- lmer(duration ~ ictus + stressed + ictus*stressed + quantity+  (1|song) + (1|performer), data = strictus_data, REML = FALSE)
summary(lindesignmod)

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: duration ~ ictus + stressed + ictus * stressed + quantity + (1 |      song) + (1 | performer)
   Data: strictus_data

     AIC      BIC   logLik deviance df.resid 
 -1779.7  -1744.7    897.9  -1795.7      581 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.1657 -0.6134 -0.0816  0.5322  3.6544 

Random effects:
 Groups    Name        Variance  Std.Dev.
 song      (Intercept) 0.0016929 0.04115 
 performer (Intercept) 0.0001313 0.01146 
 Residual              0.0026248 0.05123 
Number of obs: 589, groups:  song, 9; performer, 3

Fixed effects:
                           Estimate Std. Error t value
(Intercept)                0.221105   0.018048  12.251
ictusoff                  -0.024893   0.010000  -2.489
stressedstressed           0.034834   0.010357   3.363
quantity2                 -0.050977   0.004445 -11.469
ictusoff:stressedstressed -0.019006   0.016468  -1.154

Correlation of Fixed Effects:
            (Intr) ictsff strssd qntty2
ictusoff    -0.491                     
strssdstrss -0.467  0.896              
quantity2   -0.073 -0.073 -0.180       
ictsff:strs  0.298 -0.642 -0.663  0.169

plot(lindesignmod)

linperfmod <- lmer(duration ~ ictus + stressed + ictus*stressed  + quantity + (1|performer), data = strictus_data,REML = FALSE)
summary(linperfmod)

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: duration ~ ictus + stressed + ictus * stressed + quantity + (1 |      performer)
   Data: strictus_data

     AIC      BIC   logLik deviance df.resid 
 -1625.9  -1595.2    819.9  -1639.9      582 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.8450 -0.6461 -0.1630  0.5188  3.4975 

Random effects:
 Groups    Name        Variance  Std.Dev.
 performer (Intercept) 0.0006242 0.02498 
 Residual              0.0035523 0.05960 
Number of obs: 589, groups:  performer, 3

Fixed effects:
                           Estimate Std. Error t value
(Intercept)                0.231796   0.018150  12.771
ictusoff                  -0.044255   0.011444  -3.867
stressedstressed           0.014604   0.011858   1.232
quantity2                 -0.050451   0.005155  -9.787
ictusoff:stressedstressed  0.024063   0.018602   1.294

Correlation of Fixed Effects:
            (Intr) ictsff strssd qntty2
ictusoff    -0.567                     
strssdstrss -0.539  0.893              
quantity2   -0.086 -0.068 -0.178       
ictsff:strs  0.340 -0.628 -0.651  0.166

linsongmod <- lmer(duration ~ ictus + stressed + ictus*stressed + quantity + (1|song), data = strictus_data,REML = FALSE)
summary(linsongmod)

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: duration ~ ictus + stressed + ictus * stressed + quantity + (1 |      song)
   Data: strictus_data

     AIC      BIC   logLik deviance df.resid 
 -1781.7  -1751.1    897.9  -1795.7      582 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.1662 -0.6128 -0.0810  0.5304  3.6566 

Random effects:
 Groups   Name        Variance Std.Dev.
 song     (Intercept) 0.001824 0.04271 
 Residual             0.002625 0.05123 
Number of obs: 589, groups:  song, 9

Fixed effects:
                           Estimate Std. Error t value
(Intercept)                0.221258   0.017143  12.907
ictusoff                  -0.024862   0.010000  -2.486
stressedstressed           0.034868   0.010357   3.367
quantity2                 -0.050989   0.004445 -11.471
ictusoff:stressedstressed -0.019076   0.016469  -1.158

Correlation of Fixed Effects:
            (Intr) ictsff strssd qntty2
ictusoff    -0.517                     
strssdstrss -0.492  0.896              
quantity2   -0.077 -0.073 -0.180       
ictsff:strs  0.313 -0.642 -0.663  0.169

linnullIntmod <- lmer(duration ~ (1|song) + (1|performer), data = strictus_data, REML = FALSE)
n1 <- summary(linnullIntmod)
print(n1)

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: duration ~ (1 | song) + (1 | performer)
   Data: strictus_data

     AIC      BIC   logLik deviance df.resid 
 -1670.9  -1653.0    839.4  -1678.9      641 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.3952 -0.6996 -0.1327  0.6270  3.5442 

Random effects:
 Groups    Name        Variance  Std.Dev.
 song      (Intercept) 0.0018316 0.04280 
 performer (Intercept) 0.0001561 0.01249 
 Residual              0.0041343 0.06430 
Number of obs: 645, groups:  song, 9; performer, 3

Fixed effects:
            Estimate Std. Error t value
(Intercept)  0.20045    0.01633   12.27

compp <- anova(lindesignmod,linsongmod)
print(compp)

Data: strictus_data
Models:
linsongmod: duration ~ ictus + stressed + ictus * stressed + quantity + (1 | song)
lindesignmod: duration ~ ictus + stressed + ictus * stressed + quantity + (1 | song) + (1 | performer)
             npar     AIC     BIC logLik deviance  Chisq Df Pr(>Chisq)
linsongmod      7 -1876.7 -1845.4 945.35  -1890.7                     
lindesignmod    8 -1874.7 -1839.0 945.36  -1890.7 0.0199  1     0.8879

comps <- anova(lindesignmod,linperfmod)
print(comps)

Data: strictus_data
Models:
linperfmod: duration ~ ictus + stressed + ictus * stressed + quantity + (1 | performer)
lindesignmod: duration ~ ictus + stressed + ictus * stressed + quantity + (1 | song) + (1 | performer)
             npar     AIC     BIC logLik deviance  Chisq Df Pr(>Chisq)    
linperfmod      7 -1711.9 -1680.6 862.95  -1725.9                         
lindesignmod    8 -1874.7 -1839.0 945.36  -1890.7 164.84  1  < 2.2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

comp2 <- anova(lindesignmod,linnullIntmod)
print(comp2)

Data: strictus_data
Models:
linnullIntmod: duration ~ (1 | song) + (1 | performer)
lindesignmod: duration ~ ictus + stressed + ictus * stressed + quantity + (1 | song) + (1 | performer)
              npar     AIC   BIC logLik deviance  Chisq Df Pr(>Chisq)    
linnullIntmod    4 -1670.8 -1653 839.43  -1678.8                         
lindesignmod     8 -1874.7 -1839 945.36  -1890.7 211.88  4  < 2.2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

comp3 <- anova(linsongmod,linnullIntmod)
print(comp3)

Data: strictus_data
Models:
linnullIntmod: duration ~ (1 | song) + (1 | performer)
linsongmod: duration ~ ictus + stressed + ictus * stressed + quantity + (1 | song)
              npar     AIC     BIC logLik deviance  Chisq Df Pr(>Chisq)    
linnullIntmod    4 -1585.8 -1568.3 796.92  -1593.8                         
linsongmod       7 -1781.7 -1751.0 897.85  -1795.7 201.87  3  < 2.2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

ggplot(strictus_data,aes(ictus,euclid)) + 
  geom_jitter(color = "violet", alpha = 0.5) +
  facet_wrap("song") +
  labs(x = "beat position", y = "vowel dispersion") + 
  theme_minimal()

ggplot(strictus_data,aes(stressed,euclid)) + 
  geom_jitter(color = "orange", alpha = 0.5) +
  facet_wrap("song") +
  labs(x = "word stress", y = "vowel dispersion") + 
  theme_minimal()

ggplot(strictus_data,aes(ictus,euclid)) + 
  geom_jitter(color = "violet", alpha = 0.5) +
  facet_wrap("performer") +
  labs(x = "beat position", y = "vowel dispersion") + 
  theme_minimal()

ggplot(strictus_data,aes(stressed,euclid)) + 
  geom_jitter(color = "orange", alpha = 0.5) +
  facet_wrap("performer") +
  labs(x = "word stress", y = "vowel dispersion") + 
  theme_minimal()

NA
NA

ggplot(strictus_data,aes(x=ictus ,y=euclid,color = ictus)) + 
  geom_violin() + stat_ydensity() + theme_minimal() +
  labs( xlab = "beat" ,ylab = "euclidean distance") 

  #facet_wrap("quantity", labeller = label_both)

interaction.plot(strictus_data$stressed,strictus_data$ictus,strictus_data$euclid,xlab="word stress",ylab="vowel dispersion",trace.label = "beat")

interaction.plot(strictus_data$quantity,strictus_data$ictus,strictus_data$euclid,xlab="syllable quantity",ylab="vowel dispersion",trace.label = "beat")

ggplot(strictus_data,aes(x=stressed ,y=euclid)) + 
  geom_violin() + stat_ydensity() + theme_minimal() +
  labs( xlab = "word stress" ,ylab = "euclidean distance") 

ggplot(strictus_data,aes(x=quantity ,y=euclid)) + 
  geom_violin() + stat_ydensity() + theme_minimal() +
  labs( xlab = "syllable quantity" ,ylab = "euclidean distance")